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Lieb-Liniger model describes bosons with contact interactions in one-dimensional space. In the 
limit of weak repulsive particle interactions, there are two types of low lying excitation spectrum. 

The first is reproduced by the Bogoliubov dispersion relation, the other is believed to correspond to 
dark soliton excitations. While there are various evidences that the type II spectrum is related to 
dark solitons, it has not been shown that measurements of positions of particles reveal dark soliton 
density profiles. Here, we employ the Bethe ansatz approach and show that dark solitons emerge 
in the measurement process if the system is prepared in an eigenstate corresponding to the type II 
spectrum. We analyze single and double dark solitons as well as weak and strong interaction regime. 
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Solitons are solutions of non-linear wave equations that 
do not change their shapes in the time evolution [I] . Soli¬ 
tons are intensively studied in non-linear optics but they 
can also be created by matter waves. Ultra-cold bosonic 
atoms form a Bose-Einstein condensate (BEC) where all 
particles occupy the same single particle wave function. 
The latter is accurately described by the Gross-Pitaevskii 
equation (GPE) [0. Such a single-particle non-linear 
equation possesses bright and dark soliton solutions in 
one-dimensional (ID) space for attractive and repulsive 
particle interactions, respectively ii- Both kinds of the 
solitons have been realized experimentally 00. There 
are also attempts to create dark solitons in superfluid 
Fermi systems 110 lj|. 


Solitons constitute broken-symmetry states of many- 
body atomic systems. In the absence of an external 
trapping potential, a many-body Hamiltonian is invari¬ 
ant under spatial translation of all particles by the same 
vector. Consequently system eigenstates are also eigen¬ 
states of the corresponding spatial translation operator 
and a plot of the reduced single particle probability den¬ 
sity cannot reveal any solitonic profile because it is uni¬ 
form in space. However, even if a system is initially 
prepared in an eigenstate, solitons can emerge due to 
symmetry-breaking perturbation, e.g. measurements or 
atomic losses. This phenomenon is relatively easy to an¬ 
alyze for attractive particle interactions where the bright 
soliton formation is the result of a symmetry breaking in 
the system prepared initially in the ground state [13Ml7l| . 
On the contrary, dark solitons are related to collectively 
excited states of a repulsively interacting atomic gas. The 
question if there exist many-body eigenstates which are 
able to reveal dark solitons in the course of measurements 
of particle positions is highly non-trivial 18 -20. 


The Lieb-Liniger model, that describes bosons in ID 
space with contact interactions, can be solved by means 
of the Bethe ansatz 23 20. In the limit of weak repul¬ 


sive particle interactions, there are two branches of low 
lying excitation spectrum 0. The type I spectrum is re¬ 
produced by the Bogoliubov dispersion relation while the 
type II spectrum, related to the so-called one-hole excita¬ 
tions, is believed to correspond to dark soliton excitations 
due to the coincidence with the soliton dispersion relation 
H-El. Recently Sato et al. 011 have shown that 


superpositions of the type II eigenstates result in states 
for which the reduced single particle probability densi¬ 
ties reveal shapes very similar to dark soliton profiles. 
These states are not stationary and the soliton profiles 
diminish in the time evolution. In the present publica¬ 
tion we concentrate on single eigenstates and analyze the 
long-standing problem whether dark solitons emerge in 
the course of measurements of particles positions if the 
system is prepared in a type II eigenstate. 

We consider N bosons with repulsive contact interac¬ 
tions in a ID space with periodic boundary conditions. 
The second quantized version of the system Hamiltonian 
reads, 


H = 



d x ^d x 4> + , 


(1) 


where if is the bosonic field operator, L is the system 
size and c > 0 determines strength of the repulsive inter¬ 
actions. We use units where 2m = h = 1, with m the 
particle mass. The system is characterized by a param¬ 
eter 7 = ^ [23], where n = ^ is the particle density. 
For 7 -C 1 one deals with the weak interaction limit. 
The Bethe ansatz allows us to find eigenvalues and eigen¬ 
vectors, |{&}at), of the Hamiltonian 0 which are deter¬ 
mined by N parameters, {fc}jv = {fci,..., fcjv}, that sat¬ 


isfy the Bethe equations kjL +2 arct an (^ kj c km ^j = 

2nlj. The I 3 ’s are arbitrary integers (half integers) for 
odd (even) N which unambiguously determine an eigen¬ 
state [25|, [30 • The sets {k}w are often called the sets 
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of quasimomenta because eigenvalues of the total mo¬ 
mentum -P({fc}Ar) = X^jLi kj — 77 Ij an d the en¬ 

ergy eigenvalues -E({fc}jv) = XqLi The ground state 
of the system corresponds to the ordered set {/}at = 
{—+ 1, • ■ •, ^T}. The type II excitation 
spectrum is obtained by removing one of Ij and adding 
— 1 or + 1 [24|, HH, [3 ■ This can be called 


iV-l 

2 - 1 - 2 

one-hole excitation because one of the consecutive inte¬ 
gers (or half integers for even N ) becomes missing and 
a new integer appears at the beginning or at the end of 
the initial sequence {I}n- 

Single soliton. Let us consider even N. The ground 
state of the system corresponds to the total momentum 
P = 0. If in the sequence {L}at, we remove i and add 
;YY1 . we obtain the type II eigenstate, which we denote 
by |sol), corresponding to the total momentum P so i = 
^ 7 . The reduced single particle probability density is 
uniform, i.e. pi{x\) = jf{sol\^(xi)'^(xi)\sol} = j-. 
Thus, the result of the measurement of the position of 
a single particle is equally probable in the entire space. 
Suppose that we have performed such a measurement and 
then ask about the probability density ^ 2 (^ 2 ) for a choice 
of the position X 2 of the second particle. More generally 
we need the conditional probability density for a choice of 
the position of the j-th particle provided (j — 1 ) particles 
have been already measured (36, [371. Such a probability 
density reads 


{k} N ~j 


\ T j(xj,{k} N -j )\ S 
N-j + 1 


( 2 ) 


where the sum runs over all sets of {k}isr-j, he. over 
all eigenstates of the system with the total number of 
particles equal N — j, and 


TjOl?) j) — 'y ] Tj—1,{g}jv—j+i) 

{9l«-j+i 

x e i [ f> ({9}Ar-j+i)-f({ fe }jv-j)]3:j 

x ({k} N - j \'ip(0)\{q} N - j+1 }, (3) 


with 


Ti (a*, {</}*-!)) = 

x ({9 }jv-i|^(0)|^o), (4) 


where \tpo) is an initial eigenstate and Pq the correspond¬ 
ing eigenvalue of the total momentum. We will start with 
|^o) = | sol). The last terms in (O and (JH), so called form 
factors of the field operator, can be easily numerically 
calculated by means of the determinant formulas 38], |39| 
0 . 


Before we analyze results of simulations of the detec¬ 
tion process let us describe what we can expect from 
the mean field description. For L —> 00 , the function 
y^ntanh represents a single dark soliton solution 

of the GPE with the notch localized at x = Xq. The 



FIG. 1: (color online) (a)-(c) Results of the particle posi¬ 
tions measurements in the system prepared in the eigenstate 
corresponding to the one-hole excitation with the total mo¬ 
mentum per particle P/N = J, for N = 8 , c = 0.08 and 
L — 1, i.e. 7 = 0.01. (a) shows an example how the con¬ 
ditional probability density pj ( x ) (with values of j indicated 
in the panel) changes after successive measurements of the 
particle positions. The probability density pi (a;) for the mea¬ 
surement of the first particle is uniform as expected for sys¬ 
tem eigenstates. However, in the course of the measurements, 
pj{x) starts revealing a dark soliton notch and the last den¬ 
sity p 8 (x) becomes hardly distinguishable from the mean-field 
(MF) dark soliton profile (black solid line) corresponding to 
the same average particle momentum as P/N = f; — note 
that the mean-field profile has been shifted so that the soli¬ 
ton position xo matches the minimum of ps(x). (b) Phase of 
the mean-field solution (black solid line), whose probability 
density is shown in (a), and the phase of the wavefunction of 
the last particle (red dotted line) corresponding to the mea¬ 
surement realization presented in (a). Different simulations of 
the detection process lead to ps,{x) with the notch localized at 
different points xo■ When one shifts the obtained particle po¬ 
sitions so that xo always equals L, the average particle density 
can be prepared and it is shown in (c) (green stairs dash line) 
together with the mean-field dark soliton profile (black line). 
The average density bases on 10 4 simulations of the measure¬ 
ment process. Deviations between the average density and 
the mean-field prediction seem to be due to small N effects. 
That is, before the conditional probability density pj (x) con¬ 
verges to a final form, particles which are chosen according 
to different probability distributions distort the average den¬ 
sity. Red stairs line in (c) shows the histogram where only 
positions of last particles are collected which agrees very well 
with the mean-field soliton profile, (d)-(f) are related to sim¬ 
ilar data as (a)-(c) but for the total momentum per particle 

P/N= If- 


size of the notch is determined by the healing length 
£ = -y=- For L < 00 and £ comparable with L , dark 
soliton profile is modified. Different soliton solutions, 
that fulfill periodic boundary conditions and are char¬ 
acterized by different average values of the particle mo¬ 
mentum, are given in terms of the Jacobi functions 0 - 
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42| . There is a single-soliton solution with the average 


momentum that coincides with the momentum per par¬ 
ticle for the \sol) eigenstate, i.e. with P so i/N = J. The 
corresponding probability density is plotted in Fig. Oja) 
together with a typical example of the conditional proba¬ 
bility densities pj (x) obtained in the course of successive 
measurements of positions of particles that are prepared 
in the type II eigenstate \sol) for TV = 8 , c = 0.08 and 
L = 1, i.e. 7 = 0.01 (in this regime the mean field solu¬ 
tion is not very different from a single particle solution). 
While pi ( x ) is uniform, the successive measurements re¬ 
veal a clear notch in the probability densities. Typically, 
the measurement of the positions of 3 particles is suf¬ 
ficient to obtain Pj(x) which changes little for the rest 
of the measurements and which is similar to the mean 
held dark soliton profile; see Fig. [Ha) — note that the 
mean-held prohle has been shifted so that the soliton po¬ 
sition xo matches the minimum of pjy(x). In Fig. [TJb) 
we show plots of the phase of the mean-held solution and 
the phase of the wavefunction of the last particle cor¬ 
responding to the measurement realization presented in 
Fig. Ha). The agreement is perfect and proves that the 
detection process drives the many-body state to the dark 
soliton wavefunction. 


For N = 8 , a single simulation of the measurement 
sequence produces 8 particles positions only. In order 
to obtain an average particle density one has to perform 
many simulations. However, each simulation reveals the 
probability density pjv(x), for the choice of the last par¬ 
ticle, that is localized around a different point in the ID 
space. From a plot of pn{x) we can infer at which xo the 
anticipated dark soliton is localized. Then, all results 
can be shifted so that Xo always coincides with -§■ and 
the average particle density can be prepared as shown in 
Fig. He)- There are some deviations between the aver¬ 
age density and the mean-held dark soliton prohle which 
seem to be due to small N effects. That is, in the course 
of the measurements, a minimum of Pj(x) localizes as far 
as possible from points where the previous (j — 1 ) par¬ 
ticles have been detected. Particles chosen before Pj(x) 
converges to a final prohle, have noticeable contribution 
to the average density if N is as small as 8 . In Fig. He) 
we show also the average density of the last particles, i.e. 
those chosen according to the hnal pjy(x), which nearly 
perfectly matches the mean-held dark soliton prohle. 

For even N we have analyzed the results of the mea¬ 
surements of the system prepared in the type II eigen¬ 
state with the total momentum per particle P/N = f . 
There are other one-hole excitations corresponding to 
the momentum per particle \P\/N € {jjf, ■ ■ ■, 77 }- 
Mean-held single-soliton solutions, with the average par¬ 
ticle momentum in the same range, describe solitons 
which are rather grey than fully black because the den¬ 
sity prohles do not reach zero at soliton positions. Es¬ 
pecially the mean-held solution with the average particle 
momentum approaching 77 corresponds to the uniform 



FIG. 2: (color online) Similar results as in Fig. H a ) and 
Fig. H c ) but for the eigenstate with the total momentum 
P = 0 corresponding to two-hole excitations of the N = 8 
system with c = 0.08 and L = 1. Histograms in (b) have 
been prepared by shifting particle positions obtained in the 
simulations so that one of the two minima of pjv(x) always 
coincides with -j. 


density without any soliton signature, i.e. soliton solu¬ 
tion degenerates into a uniform superhow state |30;, 311. 


In Fig. Hd)-(f) we show the results of the particle posi¬ 
tion measurements performed in the N = 8 system that is 
prepared in the type II eigenstate with the total momen¬ 
tum per particle P/N = |j-. The mean-held approach 
predicts a grey soliton prohle which agrees very well with 
the results of the simulations, especially with the average 
density of the last particles detected in each measurement 
sequence. 


Double soliton. We may suspect that creation of a 
two-hole excitation allows us to prepare the system in an 
eigenstate which is able to reveal double dark soliton pro¬ 
hle in the course of the particle positions measurements. 
In Fig. [2] we show the results corresponding to such an 
eigenstate where in the sequence of half-integers {/}at 
that dehnes the ground state, we remove — / and ^ and 
add — | and |. Then, the total momentum of the system 
remains zero as in the ground state. The results of the 
measurements agree very well with the mean-held double 
dark soliton prohle corresponding to the vanishing aver¬ 
age particle momentum. If one analyzes the phase of the 
wavefunction of the last particle in the measurement re¬ 
alization presented in Fig. H a ) it turns out it performs 
two discontinuous jumps by 7 r at the soliton positions and 
perfectly matches the behaviour of the mean-held phase. 


Strong interaction regime. So far we have investigated 
the system in the weak interaction limit, i.e. for 7 = 0 . 01 . 
Now we switch to the strong interaction regime where the 
mean-held description is not valid. In this limit one- and 
two-hole excitations are defined in the same manner as 
we have done for the weak interactions. Examples of the 
particle positions measurements, for 7 = 1 (i.e. c = 8 , 
N = 8 ), are presented in Fig. [3] they correspond to the 
two eigenstates analyzed in the weak interaction regime. 
The conditional probability densities Pj(x) are not as 
regular as for the weak interactions. In the case of the 
one-hole excitation corresponding to the total momentum 
P = ^5, they converge to functions which show single 
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FIG. 3: (color online) Similar data as in Figs. [T][2] but for 
strong particle interactions, i.e. for c = 8, N = 8 and L = 1. 
(a)-(c) correspond to the one-hole excitation with P = -^5. 
(d)-(f) are related to the two-hole excitation with P — 0. His¬ 
tograms shown in (c) have been prepared by shifting particle 
positions obtained in the simulations so that the global min¬ 
imum of pn(x) always coincides with L. Histograms in (f) 
base on particle positions shifted so that one of the two min¬ 
ima of pn{x), as e.g. visible in (d), always coincides with -j. 
Note that black lines in (c) and (f) are not related to mean- 
field prediction as in Figs. |T][2] but to Pn(x) averaged over all 
(i.e. 10 4 ) measurement realizations. 

minimum at a certain position xq where pn{xo) ~ 0. The 
phase of the wavefunction of the last particle in the mea¬ 
surement realization presented in Fig. [3Ja) reveals jump 
at the notch position very similar as in the weak inter¬ 
action case, cf. Fig. [ljb). In the case of the two-hole 
excitation with P = 0, there are two pronounced minima 
of pjv(x). The phase of the wavefunction of the last par¬ 
ticle behaves very similarly as the phase of the double 
soliton solution analyzed in the weak interaction case, 
i.e. it performs jumps by 7 r at the minima of pn{x), 
see Fig. 03e). Distance between the two minima is no¬ 
ticeably different in different realizations of the detection 
process. Thus, contrary to the weak interaction case, the 
relative distance between the two minima of pn(x) has 
probabilistic character and cannot be considered as a well 
defined classical quantity. When we prepare the average 
particle density shifting the obtained particle positions 
so that one of the two minima of pn{x ) always coincides 
with 0.25L, only one notch is clearly visible. The other 
notch, located around x = 0.75 L, is broad and shallow. 
The density notches visible in Figs. [3tc) and (f) are nar¬ 
rower than the notches in the weak interaction limit but 
much wider than expected from the corresponding heal¬ 
ing len gth £ = Similar effect has been also observed 
in Ref. [34J where superposition of the type II eigenstates 
allows the authors to create a notch in the reduced single 
particle probability density. 


In the numerical simulations it is necessary to trun¬ 
cate sums over system eigenstates in ([2|) and @. Strong 
interaction regime constitutes the most difficult case in 
the computer simulations. In this regime, for N = 8 we 
have collected form factors ({fc} 7 |'i/S( 0 )|V>o) corresponding 
to ordered sets {/i,... ,Iy} where I\ > —9 and Iy < 9. 
Maximal values of the total momentum and quasimomen¬ 
tum of the resulting eigenstates are 5.7 kp and 6.2 kp, re¬ 
spectively, where kp is the maximal quasinromentum re¬ 
lated to the ground state of the system for N = 7. Then, 
to calculate the reduced single particle density p\(x) we 
keep only those form factors whose absolute values are 
greater than a = 10 -3 times the maximal absolute value 
of the collected form factors. All these restrictions allows 
us to obtain integrals of the reduced single particle den¬ 
sities p\(x)dx which deviate from the unity less than 
10~ 3 . In the weak interaction limit, and for a = 10 -4 , 
/ ( f pi(x)dx equals 1 within the double precision. 

In summary, we have considered the Lieb-Liniger 
model and concentrated on an analysis of the particle po¬ 
sitions measurements. Eigenstates of the system related 
to hole excitations (the so-called type II spectrum) are 
believed to correspond to dark soliton states in the weak 
interaction limit. We have shown that despite the fact 
that the corresponding reduced single particle probability 
densities are uniform, the measurements induce breaking 
of the translational symmetry and result in emergence 
of dark solitons which are localized at different positions 
in space in different realizations of the detection process. 
In the strong interaction regime, the results of the mea¬ 
surements also reveal signatures of notches in particle 
distributions. Their size is smaller than in the weak in¬ 
teraction case but greater than the healing length of the 
system. Interestingly, in the case of a two-hole excita¬ 
tion, two notches can be visible but their relative distance 
fluctuates from one realization of the detection process to 
another. 
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